A study on a self-organized criticality in a dynamical many-body system 
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A novel mechanism for the generation of self-organized criticality (SOC) is discussed in terms of 
the coupled-vibration model where the total system is forced under the uniform expansion of the 
Hubble type. This system shows a robust SOC behavior while the maximum size of the fluctuation, 
number of correlated particles in it and the temporal size of the system evolve as a function of time. 
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Observation of the ubiquitous manifestation of fractal 
structure in nature is one of the most fascinating findings 
in non- linear science 0. One direction for the under- 
standing of the origin of its manifestation was offered in 
|^|3|] with the concept of self-organized criticality (SOC). 
It was argued there that the scaling properties of the at- 
tractor is insensitive to the parameter of the model, which 
shows a robustness of the criticality. Since then, many 
models of cellular-automata were studied 0, where 
another characteristic feature is that the external driving 
of the system is much slower than the internal relaxation 
processes. 

The separation of external and internal time-scales, 
however, is not always a guiding principle for the pat- 
tern formation. Turbulent flow, for example, seems to 
have more dynamical origin. In this respect, we will seek 
a SOC state where the external driving of the system has 
a comparable time scale to that of the internal motion. 
The disappearance of the characteristic length scale and 
the robustness are used to characterize SOC. As a typical 
example which shows this feature, we propose a homo- 
geneously expanding lattice model which was motivated 
by our previous study We define a one-dimensional 
version of the previous modelQ and make an anatomy 
of SOC formation when the many-particle system is put 
in a mechanically unstable state. 

Our model is based on an expanding one-dimensional 
coupled-vibration model. The chain is composed of iVo+2 
particles of equal mass m, coupled with the nearest neigh- 
bor particles by a finite-depth potential. This chain is 
forced to obey a uniform expansion of the Hubble-type 
Physically important feature is that as the lattice 
spacing increases, each particle at the lattice point be- 
comes unstable against the motion to right or to left, 
forming a cluster. Although the following discussion will 
be applied to rather general class of interactions, we take 
Lennard- Jones (L-J) potential as a prototype 

We analyze a motion of mono-particle one-dimensional 
lattice with a lattice spacing R larger than its equilibrium 
separation uq as is given schematically in FIG.l. The 
expansion of the lattice is described with the relation. 
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FIG. 1: Sketch of the extended chain used. At equilibrium, 
random and small fluctuations are assumed. 



where the Hubble constant h is introduced. This condi- 
tion is satisfied by putting the initial velocities to each 
particle together with the corresponding velocities to the 
boundaries to which — th and (A^o + ^)th particles are 
held fixed. The equation of motion of this system is writ- 
ten as 



= V'iR + yn+i - yn) - V'{R + y„ - 2/„-i), (2) 



where yn stands for a displacement of n-th lattice par- 
ticle from its lattice point and V' is the derivative of 
the inter-particle potential. Suffix n runs from 1 to A^o 
in this equation. We assume this potential V is given 
by the Taylor polynomial expansion around R for small 
displacement r. 



V\R + r) « V'{R) + V"{R)r. 
Inserting this form to (O, we obtain 



K{2yn - yn+i - 2/«-i), 



(3) 



(4) 



R{t) = ao(l + ht) 



(1) 



where we put k for —V"{R). It should be noted here that 
the K is always negative when we deal with the stationary 
oscillation where R — ao- On the other hand, for L-J 
potential, it becomes positive when R > l.llao. In the 
following discussion, we mainly treat this situation and 
thus we assume k is positive. 



2 



To solve Q), we assume the following standard form of 
the solution, 



Un = a e:x.p{ikn =F ut). 
Inserting this form to eq.l^J, we get, 



2 4k . 2 
UJ = — sm — . 



(5) 



(6) 



This dispersion relation has the same form as for the 
normal-mode analysis for the equilibrium oscillation ex- 
cept for the fact that the motion of the lattice particle 
y„ is not oscillatory but divergent as is clear from eq.©. 
When we calculate the second derivative of in |0J, the 
following derivative appears, 



d{ujt)/dt = UJ + duj/dt X t. 



(7) 



For a non-singular form of lo, the second term in the rhs 
is smaller than the first term for sufficient small value 
of t. Therefore, we can justify the dispersion relation © 
under the condition y„ ^ R{t) and for short-time regime. 

The solution of the equation Q for small displacement 
and for small t, under the boundary condition. 



is written as 



2/0 = = VNo+i 



Dn = aexp(=Fa;i) sin(fcn), 



(8) 



(9) 



where a is a constant. In this form, the boundary condi- 
tion for n=0 is automatically satisfied and for n=A^o+li 
we obtain a relation. 



l,2,...iVo) 



No + l 

For convenience, we rewrite (llUf) and ((SJ in the form 
kii — 



(10) 



Na + 1 ' 



(11) 



^2 ^ 4« gi^2 /c. 



2- (i'=l,2,...iVo) (12) 
The general solution is written down in the form 

Un — "^-^nuQui (13) 

■^nu = \l w„\i sin JJ"^-, , (14) 

(15) 



No + l No + l ' 

a^[exp(w^i) -I- exp(-Wi,t)]/2, 



where Oi/ and are free parameters that fit the 2A'o 
initial conditions for iVo particles. 

In the v-ili normal-mode for the longitudinal motion 
(identify nR as x-coordinate), the following number of 
particles n^, are compressed or stretched as a whole, 



{No + 1) ^ K{t) 
V 2R{t) ' 



(16) 



where is a wavelength for the v-th mode expressed as 
{2L{t))lv by using the total length of the chain L{t) at 



time t. Therefore, the correlation of particle is pro- 
duced by the i/-th normal mode. The growth of the i^-th 
normal mode is governed by H15|) and is exponentially di- 
vergent. The time ordering of the growth is known from 
eas. (|ll|) and (|12|1 . For larger the ki, value becomes 
larger , which leads to larger Uv value. It means that 
the exponential growth of the fluctuation occurs step- 
wise, from the highest frequency to the lower frequency, 
in other word, from the shortest wavelength to the larger 
wavelength. When we start from the random small fluc- 
tuations of the particles' position or velocities or both, it 
corresponds to the fact that in (|15|l stands for a white 
spectrum. 

Above discussion is applicable only for small fluctua- 
tion and for short-time behavior. In other word, it is 
a linearized discussion because we assumed no coupling 
between various normal modes. For the analysis of SOC, 
however, we have to go into flnitc-time and nonlinear 
regime. We first postulate the following three hypothe- 
ses: 1) Even at finite-time, the formation of cluster starts 
from small particle number to large particle number. 2) 
At the time when Ui, correlation just starts, its spatial 
size 5c{t) is characterized by the relation 



5c{t) = \u{t)/2^7i^{t) X R(t) 



(17) 



3) The spatial size of the fluctuation for particles is 
approximately conserved even at later time. 

The hypothesis 1) and 2) are taken from the linear 
regime from (|16|l . The hypothesis 3) means that the 
size of the clusters are frozen at the time of its forma- 
tion, without stretching or shrinking at later time. It 
means that at the start of the correlation and later 
on, riy particles are trap in the mutual attractive inter- 
action bloc, free from the uniform expansion to which 
riy particles followed before. Since we assume no energy 
dissipation, the time-averaged spatial size of the bloc is 
conserved. It should be noted that when correlation 
starts, correlations of smaller particles are preserved in 
the rij, correlation, forming a nesting-boxes structure. 

Next task is to determine the dependence of the 5c{t) 
and n^it) on R{t) which satisfies H17|l . We postulate the 
following relation to generate the solution of ifTTjl . 



dSc{t) 



(18) 



where 7 stands for a non-dimensional constant. The form 
of Ihs was adopted because we study a growth process of 
rii, and 5c{t)- The form of rhs is assumed because the Ihs 
has the dimension of [length]"^ and we know that the 
system has a lattice unit parameter R{t) which has the 
dimension of [length]. The parameter 7 is a function of 
the Hubble constant h in such a way that 7-^0 for large 
value of h and 7—^1 for small value of h. The former 
condition corresponds to an enlarged copy of the system 
discussed in and the latter to a system where there 
happens one break of the chain and the left and the right 
matters move to the left and to the right without chang- 
ing their sizes. Thus the assumption of (|18|l is that Ihs 
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FIG. 2: Relation of 5c{t), N{t) and R{t) in the graph of log- 
log sale for the box size 5 and count number A'^. Cross-over 
effect and the maximum fluctuation are shown. 

is proportional to the bulk ratio cf. (|17|l . times 7, 

which stands for the amount of tension due to the com- 
petition of expanding motion and the internal attraction. 
The condition that this 7 is a constant is the condition of 
realizing the SOC, because we will see in FIG. 3 that the 
relation (|18|l is found to be well satisfied in the numerical 
calculations. 

Using (|17|) . the solution of (|18|) is written as 

ni,{t) — const X 6c{ty ■ (19) 

This rii, is directly related to the box-counting number 
N{t) for the size of (5c, as is given by the relation, 

N{t) = No/nu{t) = NaScity /const, (20) 

which is clear from (|17|) . The relation H19|l means that 
the correlation length is expressed as a power of the dis- 
tance, a typical feature of phase transition point. Here, 
n,y is considered as distance measured in unit of lattice 
spacing R(t). Since 7 depends on the expansion speed, 
the criticality is not universal here but SOC is well satis- 
fied because the system without fail goes into the critical 
state that satisfies the power- law relation (|19|) . The fact 
that all the complexity of the non-linear effect is simply 
written in a form H18|) is astonishing. 

Inserting (|19|l into (|17|l . we can get the explicit ex- 
pression of 6c{t) in terms of R{t). We can also obtain 
6c{t) in terms of the relation among 6c{t), N{t) and R{t) 
with the help of FIG. 2 drawn in log-log scale with the 
abscissa corresponding to box size S and the ordinate, to 
the count number TV. In this figure, the trace of H2(J|) 
is shown by the solid line N{t) up to a specific time t. 
For a fixed time, the maximum size of the fluctuation 
is Sc{t), beyond which the box counting number should 
follow the one-dimensional line Nref(t) shown in the fig- 
ure. That is, we define 6c{t) and N{t) = No/n,^{t) as x- 
and y-coordinates of this cross point. The line Nref{t) 



should pass through the count of unity when the box size 
6 coincides with the total chain length L, 

\ogNref{t) = -/3(log^ - log£(t)), (21) 

where /? stands for 1 in our case and 3 in the 3- 
dimensional case. 

The value of Sc{t) is the x-coordinate of the cross point 
of the line (|20|l and the line l|21|) which is expressed as 

Scitf^'* = const X R{tf, (22) 

which has the same form as the short time behaviour seen 
in dynamical scaling f^. For general spatial dimension, 
relation (|17|) is replaced by 

6,{tf = n,{t) X R{tf, (23) 

where /3=3 applies to our previous calculation 

The third hypothesis we postulated is a key to get the 
fractal growth in SOC state and the normal fractal di- 
mension calculation for a fixed system. Because of this 
hypothesis, the fractal structure formed in SOC is frozen 
out. The fact that the growth of the maximum fluctua- 
tion in our model is closely related to a cross-over effect, 
is a typical feature of our dynamical SOC state. 

Now we are ready to check our three hypotheses and 
the postulation H18|) with a help of the numerical simu- 
lation for linear chain interacting with L-J potential The 
unit of length a is the distance at which the potential 
changes sign and the energy unit e is the depth of the po- 
tential. The time unit to is composed as to — ^ (j^mjA&t. 
Before physical discussion, we first check the dependence 
of the results on the number of particles A^q. We per- 
formed the box-counting method with 100 particles and 
1,000 particles. The results hardly change when we av- 
eraged over 10 independent trials for the different initial 
states. For boundary conditions, we used a fixed bound- 
ary condition (jHl in our formalism. We checked the peri- 
odic boundary condition to simulate the infinite system 
0. We also checked the results when we include only 
the nearest neighbor interaction or allow the next neigh- 
bor one. We found that such changes altogether cause at 
most 1% change of the fractal dimension value. 

In FIG. 3, we show the box counting plot for the sys- 
tem for /i=0.1, iVo=1000 averaged over 10 events for fixed 
boundary condition and for nearest neighbor interaction. 
The value of the lattice spacing is shown in these figures 
by i?(i) in unit of the equilibrium lattice spacing. As we 
see in this figure, the fractal structure is clearly seen up 
to the cross-over point, which extends as R(t) increases. 
The linearity of the slope, which is closely related to the 
assumption of constant 7 in eq. (|18|l is well realized in this 
figure. The FIG.3(d) corresponds to iZ(i)=2.92, which 
exceeds the cut-off length 2.5 of our Lennard- Jones po- 
tential . In this situation, the formation of a cluster in 
our formalism is not physically applicable and therefore, 
the 8c value of FIG. 3(c), which is just before Rit) exceeds 
the cutoff length is shown. Actually, the behavior beyond 
this &c value starts to wave from the linear line. Another 
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FIG. 3: Numerical simulation of the box-counting for the system for h=0.1, A'^o=1000 averaged over 10 events for fixed boundary 
condition and for nearest neighbor interaction. Each panel corresponds to the time sequence of the expansion of the system. 



important tendency is the gradual decrease of the slope 
parameter as time proceeds. The change of this value is 
interpreted as due to the decrease of the stiffness param- 
eter K as R{t) increases. This causes the decrease of the 
clustering force compared with the expansion induced by 
h, which leads to a smaller fractal dimension 7. 

The same degree of good fitting with two linear lines is 
obtained for h=0.2,0.3 and 0.4. As is expected, the frac- 
tal dimension becomes smaller, the larger the h values. 
For example, the slope parameters (fractal dimensions) 
corresponding to FIG.3(c) {R{t) = 2.47) are 0.67, 0.61, 
0.56, 0.51 for h^O.l, 0.2, 0.3 and 0.4. Another important 
factor which determines the slope is the strength of the 
initial fluctuation. We adopted the initial conditions that 
A^o particles are initially at rest and are located around 
the lattice site randomly within ±p% of the original lat- 
tice spacing. This condition, when we remember ea. (|15|l . 
corresponds to = I and the strength of random vari- 
able a,y is proportional to p. The results shown in FIG. 3 
corresponds to p—5. The fractal dimension for p=5 ^ 10 
is rather stable, changes less than 3 % like three dimen- 
sional calculation in p. However, when we reduce p to 



smaller values, the fractal dimension becomes smaller, 
approaching to zero for small fluctuation limit. It is be- 
cause flj, value in p5|l becomes smaller for smaller initial 
fluctuation and the start of correlation is delayed with 
respect to time, which means the decrease of the slope of 
the fractal line. The numerical calculations show essen- 
tially the same results as in the 3-dimensional case 0, 
which suggests that this kind of SOC is common to any 
dimension. 

In conclusion, our coupled linear-chain model with 
finite-range interaction, when forced under the uniform 
expansion of the Hubble type, leads to a robust SOC 
state and shows a fractal structure and a cross-over effect. 
The slope (box-counting fractal dimension) depends on 
the expansion speed and the strength of the initial fluc- 
tuations. This SOC state in the dynamical many-particle 
system differs from the standard SOC in two respects, 1) 
relative time scales of external and internal motions, 2) 
existence of a clear cross-over effect, that has a possibil- 
ity to open a new window towards our understanding of 
fractal structure in nature. 

We thank Dr. T. Yokota for useful discussions. 
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